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Preface 


My  studies  in  computer  science  have  led  me  to  believe  that  the 
primary  duty  of  a computer  scientist  should  be  in  search  of  new  algor- 
ithms that  will  enhance  computer  efficiency.  This  thesis  evolved  in 
search  of  a new  algorithm  that  could  be  used  to  reduce  a differential 
system.  It  would  not  have  been  possible  without  the  patience  and  imag- 
ination of  Dr.  John  Jones,  Jr.,  who  pointed  me  in  the  right  direction. 
The  Air  Force  is  indeed  fortunate  to  have  him  as  a team  member. 

This  work  was  sponsored  by  the  Air  Force  Flight  Dynamics  Labora- 
tory, Wright-Patterson  AFB , Ohio.  I would  like  to  thank  Major  Eric 
Linberg,  Dr.  Dan  Repperger  of  the  Aerospace  Medical  Laboratory,  Dr. 
Robert  Craig  of  the  Air  Force  Materials  Laboratory,  and  Professors 
Constantine  Houpis  and  John  D'Azzo  of  the  Department  of  Electrical 
Engineering  for  taking  the  time  to  review  this  thesis.  The  results  of 
this  paper  should  be  of  value  to  the  United  States  Air  Force  and  others 
who  need  solutions  to  the  types  of  equations  that  I have  considered. 

I. owe  a special  debt  of  gratitude  to  my  wife  and  son  who  have  been 
as  much  a part  of  my  work  at  AFIT  as  myself. 
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Abstract 


Methods  of  reducing  linear  homogeneous  differential  matrix 
equations  of  the  form  X * A(t)X  by  kinematically  similar  transforma- 
tions are  presented  in  this  thesis.  Classical  equivalence  and  similar- 
ity relationships  between  matrices  with  elements  belonging  to  the 
polynomial  domain  with  coefficients  belonging  to  the  field  of  real 
numbers  are  obtained  which  give  necessary  and  sufficient  conditions 
for  the  existence  of  solutions  to  certain  polynomial  equations.  Block 
reduction  techniques  using  kinematically  similar  transformations  lead 
to  a discussion  of  iterative  and  non-iterative  solution  methods  for  the 
Liapunov  matrix  equation,  Liapunov  differential  matrix  equation,  Riccati 
matrix  equation,  and  the  Riccati  differential  matrix  equation. 


AN  ALGORITHMIC  APPROACH  FOR  REDUCING 
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KINEMATICALLY  SIMILAR 
DIFFERENTIAL  SYSTEMS 

I.  Introduction 


The  solution  of  differential  equations  of  the  form 


A(t)  X(t) 


(1.1) 


where  A(t)  is  an  n x n matrix,  has  a wide  variety  of  applications  in 
the  engineering  sciences.  In  a computer  environment  solutions  to 
Eq  (1.1)  can  sometimes  involve  a time  consuming  iterative  technique 
with  questionable  accuracy.  Classical  methods  of  reducing  Eq  (1.1) 
to  a simpler  system  have  generally  been  limited  to  the  case  where  the 
matrix  A is  a constant  matrix  or  has  periodic  or  almost  periodic 
elements.  For  the  case  where  the  matrix  A has  constant  elements, 
reduction  techniques  have  involved  the  computation  of  eigenvectors. 
These  eigenvectors  are  sometimes  difficult  to  compute  when  the  matrix 
A has  multiple  characteristic  roots  and  the  process  cannot  always  be 
accomplished  in  an  efficient  way  on  a computer.  This  thesis  will 
develop  a new  algorithmic  approach  for  reducing  differential  systems 
with  matrices  such  as  A in  Eq  (1.1)  for  both  the  constant  and  some 
variable  cases  so  that  the  computation  of  the  eigenvectors  or  a long 
iterative  process  is  not  required  in  order  to  obtain  a solution.  The 
reduction  technique  will  allow  one  to  reduce  the  matrix  A by  blocks 
either  partially  or  completely  under  rather  general  conditions  on  A. 
This  reduction  of  Eq  (1.1)  to  a simper  but  similar  system  is  extremely 
advantageous  in  terms  of  the  computer  time  and  memory  requirements 
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needed  Co  solve  large  differential  systems. 

The  Idea  of  reducing  a matrix  to  a simpler  but  similar  form  was 
presented  by  H.  Loevy  (Ref  1).  He  defined  a relationship  between  two 
matrices  whose  elements  are  functions  of  a variable  x as  follows: 

Definition  1.1.  Let 

P - (P^ (x)),  P°  - (d  PdU  <x)),  (i,j  - 1,2,. ...n)  (1.2) 

then  we  may  write  A B (A  is  similar  to  B in  the  sense  of  Loewy)  if 
there  exists  a matrix  P whose  elements  are  differentiable  functions  of 
x and  whose  determinant  does  not  vanish  identically  such  that 

B(x)  - -PDP-1(x)  + P(x)  A(x)  P_1(x)  (1.3) 

Further  research  in  this  area  has  been  done  by  Coppel  (Ref  2) , and 
as  presented  by  Berkey  (Ref  3)  it  states  that  a differencial  system 
in  the  form  of  Eq  (1.1)  is  kinematically  similar  to  another  equation 

•Y(t)  = B(t)  Y(t)  (1.4) 

of  the  same  form,  if  there  exists  a continuously  differentiable 
invertible  matrix  function  P(t)  satisfying  the  equation 

P(t)  = A(t)  P(t)  - P(t)  B(t)  (1.5) 

which  is  bounded  together  with  its  inverse  on  the  real  line.  The 
vector  change  of  variable 

XU)  - P(t)  YU)  (1.6) 

transforms  Eq  (1.1)  into  (1.4)  so  that  they  both  have  the  same  stability 
and  boundedness  properties. 
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In  the  process  of  developing  a P matrix,  as  suggested  by  Loevy  and 
Coppel,  solutions  of  the  Liapunov  equation  and  the  Riccati  equation 
will  be  investigated.  These  equations  by  themselves  have  a wide  variety 
of  applications  in  the  field  of  engineering  sciences.  Computational 
methods  used  to  solve  these  equations  can  be  found  throughout  the 
literature.  Of  particular  note  is  work  done  by  Leuthauser  (Ref  4) 
which  will  be  discussed  in  Chapter  IV. 

In  this  thesis  capital  letters  denote  matrices  and  matrix  func- 
tions, while  lower  case  letters  represent  scalars  and  scalar  functions. 

T is  the  transpose  symbol,  and  I is  the  identity  matrix. 
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II.  Solutions  to  Polynomial  Equations 

The  main  purpose  of  this  chapter  is  to  show  how  the  solutions 

to  certain  matrix  equations  with  elements  belonging  to  a polynomial 

domain  F[z],  where  F is  the  field  of  real  numbers,  implies  the 

equivalency  and  similarity  of  pairs  of  matrices.  The  solutions  of 

the  matrix  equations  are  in  themselves  a worthy  topic  and  can  be 

found  in  various  works  such  as  Lhose  of  B.  Porter  (Ref  5),  and  as 

will  be  shown  later  in  this  chapter  as  a means  to  evaluate  the 

At 

fundamental  matrix  e where  A is  the  n x n matrix  of  the  differential 
system  given  by  Eq  (1.1).  This  in  turn  will  have  a direct  applica- 
tion to  obtaining  a fundamental  solution  to 

XU)  - a XU)  (2.1) 

where  A has  elements  belonging  to  the  field  of  real  numbers. 

2.1  Solutions  to  Certain  Polynomial  Equations 

Definition  2.1.  The  matrix  A is  said  to  be  similar  to  the 
matrix  B if  and  only  if  there  exists  a non-singular  matrix  P such 
that  P * A P = B. 

Definition  2.2.  The  matrix  A is  equivalent  to  the  matrix  B 
if  and  only  if  there  exists  non-singular  matrices  P and  Q such  that 
P A Q - B. 

Theorem  2.1.  Let 

A(z)  X(z)  = C(z)  (2.2) 

where  A,  X,  and  C are  n x n matrices  having  elements  belonging  to 
F[z],  then  the  following  2n  x 2n  matrices  are  equivalent: 
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(2.3) 


Proof.  By  definition  2.2  one  must  demonstrate  that  there 
exists  a pair  of  matrices  P and  Q so  that  the  pair  of  matrices  given 


by  (2.3)  are  equivalent.  Let  X(z)  be  a solution  of  Eq  (2.2)  and  let 


P - 


-A(z)X(z) 

1 J • 


(2.4) 


Using  these  matrices  as  required  in  definition  2.2  yields 


’i 

-A(z)X(z) 

A (z) 

C(z) 

J 

0 

r 

i J 

_ 0 

I 

Lrl 

H 

[a(«) 

C(z)  - 

A(z)X(z)] 

A(z) 

H 

L 0 

i J 

W 

0 

N 

since  Eq  (2.2)  has  a solution  by  hypothesis. 

Theorem  2.2.  If  X(z)  is  a solution  of  Eq  (2.2)  then  the 
following  pair  of  2n  x 2n  matrices  are  similar: 


A(z) 

C (z) 

A(z) 

0 

_ 0 

0 

f 

0 

o_ 

^ • By  definition  2.1  one  must  produce  a non- singular 
matrix  P so  that  the  pair  of  2n  x 2n  matrices  in  (2.6)  are  similar. 
Let 


i 

-X(z) 

p-1  = 

i 

X(z) 

0 

I 

t 

_0 

I 

(2.7) 


5 


and  using  these  matrices  as  required  in  definition  2.1  yields 


I 

X(z) 

A(z)  C(z)"' 

P 

-X(z) 

0 

[ 

0 0 

L° 

I 

|"a(z) 

C(z) 

- A(z)X(z)l 
1 

. f 

'a(z)  I 

_fl 

(2.8) 


Theorem  2.3.  Let  X(z)  be  a solution  of  the  matrix  Eq  (2.2) 
then  the  following  equation  must  hold  true: 


X(z$) 

A(z) 

C(z) 

0 

0 

-X(z) 


- (0) 


(2.9) 


Proof.  Multiplication  of  the  left-hand  side  of  Eq  (2.9)  will 
yield  the  zero  matrix  since  Eq  (2.2)  is  assumed  to  have  a solution 
X(z). 

Theorem  2.4.  Let  X(z)  be  a solution  to  Eq  (2.2)  then 


&(z)  | C(z£) 


-X(z) 


= (0) 


(2.10) 


Proof.  Multiplication  of  the  left-hand  side  of  Eq  (2.10)  will 
yield  the  zero  matrix  since  Eq  (2.2)  is  assumed  to  have  a solution  X(z). 

Before  proceeding  with  theorem  2.5  it  will  be  necessary  to  standard- 
ize some  notation: 

1)  For  polynomial  elements  of  a matrix, a^ (z) | a^ (z) , where 
z is  a real  variable,  means  that  the  element  a^(z)  of  matrix  A 
divides  into  the  element  a„  (z)  and  the  quotient  is  also  a polynomial 
in  z. 
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2)  (a^(z),  ajj(z))  means  the  lushest  polynomial  common 

factor  of  aii(z)  and  a4l(z).  For  example, 


jj 


( (z+1),  (z+1)*  ) - (z+1) 


(2.11) 


and  also 

( (z+1),  z)  - 1 (2.12) 

Eq  (2.12)  implies  that  a^(z)  “ (2+l)  is  relatively  prime  to 
ajj(z)  “ z.  A useful  result  that  will  be  used  in  conjunction  with 
this  notation  is  that  if  f(z),  g(z)  are  polynomials  in  z and  f(z), 
g(z)  are  relatively  prime  to  each  other,  that  is 

(f (z),g(z))  = 1 (2.13) 

then  there  exist  a(z)  and  b(z)  which  are  polynomials  in  z such  that 

a(z)f (z)  + b(z)g(z)  = 1 (2.14) 

Theorem  2.5.  If  the  pair  of  2n  x 2n  matrices 


A(z) 

C(z) 

A(z) 

0 

_ 0 

A(z)_ 

» 

0 

A(z) 

(2.15) 


are  equivalent  then  Eq  (2.2)  has  a solution  x(z)  with  polynomial 
elements  in  z. 

Proof.  Using  the  results  of  W.  E.  Roth  (Ref  6)  the  equivalency 
of  the  pair  of  matrices  of  (2.15)  implies  the  existence  of  non-singular 
matrices  P and  Q as  required  by  definition  2.2  such  that 


P(z)A(z)Q(z)  = diag  (a^U),  a^2(z) » • • • a^a(z) , 0,0, . . .0}  (2.16) 
where  a^i(z)|a'  (z)  for  i<j=2, 3, . . .a.  Using  the  same  matrices  P(z), 
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Q(z),  and  C(z)  yields 


P(z)C(z)Q(z)  - C' (z) 


where  C' (z)  has  elements  c^(z)  which  are  polynomials  in  z such 


c^(z)  - k1j  (z)  (a^Xz) , a^(z)),  (i, j-1,2, . . .n) 


where  k^(z)  is  a constant  or  a polynomial  in  z. 

Since  (a'  (z),  a'  (z))  = a'  (z),  Eq  (2.18)  reduces  to 
il  act  11 


c|j(z)  = k±j(z)  a^z),  (i,j=l,2,...n) 


Thus  an  element  by  element  solution  to  the  equations 


aii(z)uij(z)  = Cij^’  (i*  J=1»2»  • • •") 


would  have  a solution 

u^z)  = kij  (z) , (i,j=l,2,...n) 

In  matrix  notation  the  following  equation  would  be  valid  from 
Eq  (2.20): 

A'U  = C' 

But  from  Eqs  (2.16)  and  (2.17)  this  is 
[PAQ]U  = PCQ 

Multiplying  Eq  (2.23)  on  the  left  by  P * and  on  the  right  by  Q 
yields 


A[QUQ_1]  = C 


(2.17) 

that 

(2.18) 

(2.19) 

(2.20) 

(2.21) 

(2.22) 

(2.23) 

(2.24) 
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which  implies 


X(z)  - Q U Q 


(2.25) 


is  a solution  to  Eq  (2.2). 

Theorem  2.5  implies  that  by  finding  non-singular  matrices  P and  Q 
that  produce  the  diagonal  matrix  of  Eq  (2.16),  one  can  then  generate 
a matrix  U by  solving  Eq  (2.20)  element-wise  and  thus  produce  a solu- 
tion X(z)  to  Eq  (2.2)  by  computing  Q U Q-1  in  Eq  (2.25). 

Example  2.1.  Find  X(z)  in  Eq  (2.2)  when 


1 0 
0 z2H 

0 -z 


A(z)  = 0 z +z  0 

2 


C(z)  * z^-z  z2+z  z^-z  (2.26) 

2,  3 2 2 

z +z  z -z  -z  z 


10  0 
Oil 
0 10, 


10  0 
0 0 1 
0 1-1 


Eq  (2.16)  yields 


10  0 
A'  = 0 z 0 

0 0 z2+z 


P A(z)  Q 


Eq  (2.17)  yields 


2 

, 2 

z 

z 

l-2z 

3.  2 

3.  2 

2 

z +z 

z +z  -z 

-z  +z 

N 

CO 

1 

N 

N 

CO 

1 

N 

3 2 

-z  +z  +".z 

(2.27) 


(2.28) 


(2.29) 


Now  constructing  u^(z)  from  Eq  (2.20)  yields 


11  * 2 

2 

u12  ’ 2 

„13-l-2,2 

2 

21  ‘ 2 +2 

2 , 
u22  » z +Z-1 

“23  ' '2+1 

► 

(2.30) 

31  ‘ z‘l 

u32  _ Z-1 

u33  ’ -2+2 

Eq  (2.25)  will  yield 
z 

X(z)  = |z-l 


. 2 2 
1-z  z 


z— 1 


z2+l  z2-l  z2  J 


QUQ 


-1 


(2.31) 


It  should  be  noted  that  when  Eq  (2.2)  has  rational  polynomials  as 
elements,  it  can  be  reduced  to  the  case  where  the  elements  are  poly- 
nomials by  multiplying  Eq  (2.2)  by  the  least  common  multiple  of  the 
denominators  of  A(z)  and  C(z). 

Example  2.2.  If  Eq  (2.2)  were  of  the  form 


z2+z 

z 

z 

2* 

z 

Z2-1 

z+1 

X(z)  = 

z-1 

z+1 

z+5 

1 2 

2z 

_z2+2z+1 

z 

(2+1)  2 

z-1 

(2.32) 


then  it  could  be  reduced  to  a form  of  only  polynomial  elements  by 

2 

multiplying  A(z)  and  C(z)  by  (z-l)(z+l)  which  is  the  least  common 
multiple  of  the  denominator  terms  of  the  elements. 

The  next  theorem  will  extend  the  ideas  presented  to  the  solution 
of  a differential  equation  whose  coefficient  matrix  has  elements 
belonging  to  the  polynomial  domain  F[z]. 
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Theorem  2.6.  If  X(z)  is  a solution  of  the  differential  equation 


X(z)  =■  A(z)  X(z)  - C(z) 
then  the  following  matrices  are  equivalent 


and  thus 


(2.33) 


(2.34) 


(2.35) 


(2.36) 


The  next  theorems  use  an  approach  similar  to  the  previous  theorems 
in  order  to  solve  a matrix  equation  of  the  form 


A(z)X(z)  - X(z) A(z)  = C(z) 


(2.37) 


where  A(z),  C(z),  and  X(z)  are  n x n matrices  having  elements  belonging 
to  the  polynomial  domain  F[z], 

Theorem  2.7.  If  X(z)  is  a solution  of  Eq  (2.37)  then  the 
pair  of  matrices  given  in  (2.15)  are  similar. 


11 


Proof.  Let 


I -X(z) 


I X(z) 


p 

o 1 

I 

I ■ 

P * 

b 

! i 

in  order  t 

.o  satisfy  definition  2.1. 

Thus 

I 

X(z) 

A(z) 

C(z) 

I 

-X(z) 

0 

I 

0 

A(z) 

0 

I 

(2.38) 


(2.39) 


provided  X(z)  is  a solution  to  Eq  (2.37). 

Theorem  2.8.  If  X(z)  is  a solution  of  the  differential  equation 


X (z)  = A(z)X(z)  - X(z)A(z)  - C(z) 
then  the  following  pair  of  matrices  are  similar 


(2.40) 


A(z)  C(z)+X  (z) 


A(z) 

0 

0 

A(z) 

(2.41) 


Proof.  Let 


-X(z) 

I , 


I X(z) 


(2.42) 


in  order  to  satisfy  definition  2.1.  Thus 


I X(z)  A(z)  C(z)+X  (z)  I -X(z) 


I 

-X(z) 

0 

I 

(2.43) 

with  X(z)  a solution  of  Eq  (2.40). 

The  next  theorem  is  given  in  Nerring  (Ref  7:115)  and  will  be 
needed  to  prove  theorem  2.10. 

Theorem  2.9.  If  the  (trace  of  C(z))  = 0 then 

(trace[PC(z)P_1] ) = 0 (2.44) 

Theorem  2.10.  If  the  pair  of  matrices  in  (2.15)  are  equivalent 
then  Eq  (2.37)  has  a solution  X(z)  provided  the  following  conditions 
hold: 

1)  Trace  [C(z)J  = 0 

- 2)  There  exists  a non-singular  matrix  P such  that  P ^AP 

equals  a diagonal  matrix 

3)  For  all  o and  6,  a t 8 the  invariant  factors  a (z), 

aa 

a..(z)  of  P ^A(z)P  must  have  the  property  that  a (z)  ^ 
pp  aa 

a66(z)‘ 

Proof.  Let  P A(z)  P-1  = diag  {a;L1(z),  a22(z),  •••  a6g(z), 

0,  0,  ...  0}  where  a^z)  |a_  (z) , (j=2,...B)  and  a±i(z)  * a^  ^ (z)  to 
satisfy  condition  3.  Since  Trace  [C(z)]  = 0 then  Trace  [PC(z)P  ^]  = 0 
by  theorem  2.9.  Using  the  results  of  Roth  (Ref  6)  on  equivalent 
matrices  the  following  matrix  equations  can  be  formed: 

[ 


A(z) 

0 

0 

A(z) 

the  elements  c'^  are  then  the  highest  common  factors  of  a^(z)  and 

aj j (Z> ’ thuS 

cj^(z)  = kij(z)(a^i(z),  ajj(z)),  (i,  j-1,2, . . .n)  (2.47) 

where  k^(z)  is  a constant  or  a polynomial  in  z.  Therefore,  to  solve 
the  equations 


aoa(z)  ua6(Z>  • Ua8(z)a6B(z)  " ca6(z) 


(2.48) 


where  u^are  unknown  polynomial  functions,  one  must  investigate  three 
possibilities: 

Case  1.  o < 8 

Case  2.  a > 8 

Case  3.  a = 8 

Case  1.  For  a < 8 then  a^(z)  |a^g(z)  and  Eq  (2.48)  can  be  written 

as 

aao(z)ua6(z)  “ Ua6(z)[1a6(z)aota(z)1  “ kaB(z)aaa(z)  (2'49) 

a’(z) 

where  1 . = , . . Therefore 

a8  a'  (z) 
aa 


ua8(z)  ' Ua8(z)  1a6(z)  = ka8(z) 
ka8(z) 


“a6(z)  " 1-1ae(z) 

which  is  a ratio  of  polynomials. 


(2.50) 

(2.51) 


Case  2.  For  a > 8 then  a„~la  (z)  and  Eq  (2.48)  can  be  written  as 

” PP  OtOl 


where  1 „(z)  = 
a8 


‘8B(Z)1 

UaB(z) 

t (z) 

aa 

Thus 

!66(Z) 

aB  68 


aB  66 
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(2.53) 


kQ(3<z> 


l„p(z>  ~ 1 


(2.54) 


Case  3.  For  a = 6 then  c „ * 0 by  condition  1 and  the  solution 

" ■ ■ 1 Qp 

u a is  any  arbitrary  polynomial  in  z. 
ap 

The  elements  u . which  form  a matrix  U can  be  used  to  solve 
ap 

Eq  (2.37)  since 

A'U  - UA'  - C' 

or  in  other  words 

(PA(z)P-1)U  - U(PA(z)P-1)  = PC(z)P_1  (2.55) 

Multiplying  Eq  (2.55)  on  the  left  by  P ^ and  on  the  right  by  P yields 

(2.56) 

Therefore  a solution  to  Eq  (2.37)  would  be  given  by 


A(z)(P_1UP)  - (P_1UP)A(z)  = C(z) 


X(z)  = P-1UP 


Example  2.3.  Find  X(z)  in  Eq  (2.37)  with 


(2.57) 


1 

o 

o 

rn  4 2 , 2 3 1 

0 z -z  -z+1  z -z 

A(z)  = 

2 

0 z +z  0 

C(z)  = 

4 3 , 2 3 3 

z +z  +z  z z 

. 2 

4 . 2 3 

0 -z  z 

9 

_Z-Z  -1  2 Z 

1 0 
0 1 
0 1 


,-l 


10  0 
0 0 1 
0 1-1 


(2.58) 


(2.59) 
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w 


Forming  the  matrix  M of  Eq  (2.46)  yields 


10  0 


0 z 0 


0 0 z+z 


0 0 0 


z3-z2+z-l 

4.3  2 

z +z  -z 


z2-z3  z4+z3-2z2-z+1 

0 z3+z2 


(2.60) 


0 0 0 


0 0 0 


Solving  u 0 from  Eq  (2.48)  yields 
op 


2 2 

= arbitrary;  u^2  = z ; ■ 1 - z 


u^  * z +1;  u22  = arbitrary;  = -1  - z 


u^  = z ; u^2  = z;  = arbitrary 


(2.61) 


Solving  Eq  (2.57)  for  X(z)  yields 


0 1 z 

X(z)  = z2  z z = P-1  U P 


(2.62) 


1 l-2z  -z. 


It  should  be  noted  that  in  the  examples  generated  uag  are  polynomials 
in  z and  not  the  more  restrictive  case  implied  by  Eqs  (2.51)  and  (2.53). 
It  is  possible  that  such  is  always  the  case,  but  a proof  could  not  be 
generated  to  substantiate  this. 


The  next  example  will  show  a direct  application  of  W.  E.  Roth 


(Ref  6)  which  uses  a similar  technique  to  solve 


AX  - XA  = C 


where  A and  C belong  to  the  field  of  real  numbers. 


(2.63) 
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Examp le  2.4. 
1 0 

A - 

1 0 


-2  2 

-1  2 


(2.64) 


From  Roth’s  lemma  (Ref  6)  if  the  matrices  of  Eq  (2.15)  are  similar 
then  the  following  matrices  are  equivalent  with  polynomial  elements  in 
x with  real  coefficients: 


(A-xI) 

C 

0 

(A-xI) 

(A-xI) 

0 


0 

(A-xI) 


(2.65) 


According  to  Roth  this  then  implies  the  existence  of  matrices  X and  Y 
such  that 


(A-xI)X  - Y(A-xI)  = C 

Finding  matrices  P and  Q for  (A-xl)  as  in  theorem  2.5  yields 

0 1 lx 

p * Q 

1 x-lj  , 0 1 

Thus  for  the  example 

P(A-xI)Q  =10  = A' 

2 

0 x-x 


(2.66) 


(2.67) 


P C Q 


-1  2-y 

. 2 
-x-1  -x 


(2.68) 


(2.69) 


Forming  the  matrix  M in  Eq  (2.46)  yields 


0 x-x" 
0 0 
0 0 


2 

-x  +x 


(2.70) 
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Then  as  in  theorem  2.10  solving  an  equation  of  the  form 


A'U  - VA*  - C' 


will  lead  to  a solution  of  Eq  (2.66)  with 


(2.71) 


X(x)  *QUQ 


Y(x)  - P-1VP 


-1 


(2.72) 

(2.73) 


where 


X « X + X,x+.. .X  xh 
o 1 p 


Y ■ Y + Y,x+. ..Y  x^ 
o 1 q 


(2.74) 

(2.75) 


such  that  X^  and  Y^  are  n x n matrices  with  elements  belonging  to  the 
field  F.  Thus  using  Eq  (2.71)  to  solve  for  U and  V element  by  element 
for  the  example  yields 


(2.76) 


U11 

• #> 

o 

II 

vii  = 1; 

U12  = ^2_x^; 

v12  = 0 

U21 

= 0; 

V21  = °; 

U22  = 1; 

v22  = 0 

Therefore 

To 

-1 

[i  -1 

(2.77) 


Applying  Eqs  (2.72)  and  (2.73)  yields 
0 2 
0 1| 


X(x)  = 


Y(x) 


(2.79) 


With  Eq  (2.79)  as  a solution  to  Eq  (2.66),  a solution  to  Eq  (2.63)  is 
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X - Z AY 
i-0 


(2.80) 


or 


Z X B 
i-0 


(2.81) 


where  A°  = B°  = I 


Thus  using  Eqs  (2.80)  or  (2.81)  will  yield 


0 2 
0 1 


(2.82) 


as  a solution  to  Eq  (2.63). 

Example  2.4  shows  that  two  solutions  may  be  readily  computed  for 
Eq  (2.63).  The  next  theorem  could  prove  useful  if  two  such  solutions 
are  available. 

Theorem  2.11.  If  X^  and  X^  are  solutions  to  Eq  (2.63)  then 
X^  - X2- is  a solution  to 


AW  = WA 


(2.83) 


Proof.  Let 


AXX  - XjA  = C 


(2.84) 


also 


ax2  - x2a  = C 


(2.85) 


Subtracting  Eq  (2.85)  from  (2.84)  yields 
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(2.86) 


A(XX  - X2)  - (X1  - X2)A  - 0 


By  leCCing 

W » Xx  - X2  (2.87) 

a solution  to  Eq  (2.83)  is  found. 

2.2  A Fundamental  Solution  to  Eq  (2.1) 

This  section  will  show  how  the  solution  of  polynomial  matrix 

equations  can  be  used  to  compute  a solution  to  Eq  (2.1). 

Theorem  2.12  (Bronson  (Ref  8)).  The  matrix  system  in  Eq  (2.1) 

with  X(t  ) * k.  has  the  solution 
o 

I(t)  = {eX(t"to)]T  (2.88) 

Theorem  2.13  (Bronson  (Ref  8)).  If  A is  a matrix  having  n rows 
and  n columns  then 

eAt  = a (t)An  ^tn  L+a  (t)An  2tn  2+. . .+a  (t)At+a  ( t ) I (2.89) 

n-1  n-^  i o 


\ 


where  a (t),  a,  (t),...a  , (t)  are  polynomials  in  t which  must  be 

o 1 n-i 

determined  for  a given  A in  Eq  (2.1). 

Theorem  2,14  (Bronson  (Ref  8)).  Let  A be  as  in  theorem  2.13 
and  define 

y(X)  = a ,Xn  1 + a „Xn  2 + . . . + a, X + a (2.90) 

n-1  n-2  1 o 

then  if  X^  is  an  eigenvalue  of  At 
Y± 

e 1 = y(X±)  (2.91) 


.... 
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Furthermore,  if  A^  is  an  eigenvalue  of  multiplicity  k,  k>l  then  the 
following  equations  are  also  valid: 

Xi  _ 

dA  , 


Ai  . d2y(A) | 


dk-1  Y(A) 


(2.92) 


Example  2,5.  Find  a solution  to 

x(t)  = 1 1(0 

4 3 


(2.93) 


with  X(to)  = |_J 

At 

Applying  theorems  (2.13)  and  (2.14)  in  order  to  find  e produces  the 
characteristic  roots  of  At  as:  -t,  5t 
By  applying  Eq  (2.91)  yields: 

v(-t)  = a, (t)  ’ (-t)  + a (t)  = e-t  (2.94) 


(2.95) 


Y (~t)  = a1(t) 

(-t) 

+ aQ(t)  = e 

y(5t)  = aL(t)  • 

(5t) 

+ aQ(t)  = e 

Transforming  Eqs  (2.94) 

and 

(2.95)  into 

r ii  pi 

(t)l 

- n 

1st  iJ  L,0 

(t)J 

U5tJ 

(2.96) 
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Eq  (2.96)  is  now  in  the  form  of  Eq  (2.2)  and  theorem  2.5  can  be  used 


to  solve  for  a and  a with 
o 1 


5 1 

1 -1 

and  applying  Eq  (2.16)  yields 
6 0 
0 -6t 


0 1 
1 0 


(2.97) 


(2.98) 


Similarly  applying  Eq  (2.17)  by  expanding  the  matrix  C(t)  of  Eq  (2.96) 


to  square  dimensions  yields 
C’  = 


/C  . 5tN  -t.  5tT 

(5e  +e  ) (5e  +e  ) 


, “t  5t.  , -t  5tv 

(e  -e  ) (e  -e  ) 


Constructing  the  matrix  U from  Eq  (2.20)  yields 


(2.99) 


rc  -t  5t 

5e  -e 


U = 


, -t  5t 
5e  -e 


+6 

-t  5t 
e -e 


-6t 


+6 

e_t-e5t 
^6t~  J 


(2.100) 


then  solving  Eq  (2.25)  implies 


“x(t) 

’ 5t  -t 
e -e 

6t 

“o(t). 

r “t , 5t 
5e  +e 

. 6 J 

(2.101) 


Thus  substituting  and  into  Eq  (2.89)  yields 


At  1 


2e5t+4e_t 

4e5t-4e-t 


2e5t-2e-t 

4e5t+2e-t 


(2.102) 


Then  using  theorem  2.12,  a solution  to  Eq  (2.93)  is 
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1 

6 


(2.103) 


X(t) 


4e5t+2e 

8e5t-2e 


By  taking  advantage  of  equivalency  and  similarity  relationships,  this 
chapter  developed  methods  to  solve  equations  with  polynomials  as  elements. 
Other  methods  are  of  course  available,  in  particular  the  works  of  M.T. 
McClellan  should  be  referenced.  (Ref  9),  (Ref  10). 
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III.  Reducibility  Using  Solutions  of  the  Liapunov  Equation 


This  chapter  will  apply  the  ideas  of  kinematic  similarity  presented  in 


Chapter  I to  specific  forms  of 


X(t)  = A(t)  X(t) 


(1.1) 


The  first  section  will  show  how  reducibility  of  a differential  system 
to  a kinematically  similar  differential  system  is  applicable  and  includes 


the  classical  reduction  method.  The  second  section  will  apply  the 
reduction  techniques  to  Eq  (1.1)  for  the  constant  case  by  solving  the 
Liapunov  Equation.  The  third  section  will  extend  the  ideas  to  the  time 
dependent  case,  and  the  final  section  of  this  chapter  will  discuss 
possible  methods  that  could  be  used  to  solve  the  Liapunov  Equation. 

3.1  Reducibility  of  Eq  (1.1)  to  Upper  Triangular  Form 

Theorem  3.1  (Noble  (Ref  11)).  If  A is  a n x n constant  matrix 
with  eigenvalues  (which  need  not  be  distinct)  then  there  exists 

a non-singular  matrix  P such  that 


— 1 

P A P = U 


(3.1) 


where  U is  an  upper  triangular  matrix  with  A,...  A along  the  diagonal. 

i n 

Theorem  3.2,  If  A satisfies  the  conditions  of  theorem  3.1  then 
the  differential  system  of  Eq  (1.1)  is  kinematically  similar  to 


Y(t)  = U Y(t) 


Proof.  Let 


*(t)  = P Y(t) 


(3.2) 


(3.3) 


where  P is  as  defined  in  theorem  3.1. 


I 

L 


Differentiating  Eq  (3.3)  and  the  substitution  into  Eq  (1.1)  yields 


X(t)  - A P Y(t)  - P Y (t) 


(3.5) 


,-l 


Since  P exists,  Eq  (3.5)  transforms  to 


Y(t)  - [P-1  A P]  Y(t) 


(3.6) 


which  by  theorem  3.1  is  the  desired  result. 
Example  3.1.  Reduce 
5 4 3l 

|-1  0 -3 

1-2  1 


X(t) 


X(t) 


(3.7) 


The  matrix  P of  theorem  3.1  from  Noble  (Ref  11)  is 


10  0 

-110 
-10  1 


„-l 


10  0 
110 
10  1 


(3.8) 


Using  Eq  (3.6),  the  system  of  Eq  (3.7)  reduces  to 


Y(t) 


-2  3 4 

0 4 2 

0 0 4 


Y(t) 


(3.9) 


3.2  Reducibility  of  Eq  (1.1)  Using  a Solution  of  the  Liapunov  Equation 
Theorem  3.3.  The  differential  system  of  Eq  (1.1)  in  the  form 


X(t) 


[■  -] 


X(t) 


(3.10) 


where  A,  B,  and  C are  n x n matrices  with  elements  belonging  to  the  field 
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of  real  numbers,  is  kinematically  similar  to 


(3.11) 


by  the  transformation 


(3.12) 


where  X is  a solution  to  the  equation 


AX  - XB  = C 


(3.13) 


Proof.  Let 


X(t)  = P Y(t) 


(3.3) 


where 


0 I IJ  , 


■i  i-x- 


.0  I I. 


(3.14) 


Differentiating  Eq  (3.3)  yields 


(3.15) 


Substituting  Eq  (3.3)  into  Eq  (3.10)  and  equating  with  Eq  (3.15)  yields 


-a  i -cl  ri  I x- 


.0  I bJ  Lo  I i. 


•I  XT  t— >• 

Y(t) 

.0  I. 


(3  16) 


Multiplying  Eq  (3.16)  on  the  left  by  P yields 
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ffj 

tl 

‘A 

AX 

- C - 

XBl 

.0 

B 

J 

which  reduces  to 

Eq 

(3.11)  with  X 

Example 

3.2. 

Reduce 

-i 

1 

0 

-1 

-3 

i(t)  * 

1 

0 

-1 

-2 

0 

0 

0 

-1 

0 

0 

0 

-1 

with 

- [i 

a 

9 

B = 

R 

a solution  to  Eq 

(3.13) 

is 

Y(t) 


X(t) 


Ei  -0 


x 


c a 


thus  using  the  matrices  of  (3.14) 


1 0 
0 1 


0 0 
0 0 


1 1 
1 0 


1 0 
0 1 


1 0 
0 1 


0 0 
0 0 


-1  -1 

-1  0 


1 0 
0 1 


the  substitution  of  X in  Eq  (3.3)  reduces  Eq  (3.18)  to 


Y(t) 


0 

0 


0 0 
0 0 


0 0 
0 0 


0 -1 

0 -1 


Y(t) 


(3.17) 


(3.18) 


(3.19) 


(3.20) 


(3.21) 


(3.22) 


A slightly  different  version  of  theorem  3.3  that  will  be  used  in  Chapter  IV 


1 
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Is  as  follows: 


Theorem  3.4.  The  differential  system  of  Eq  (1.1)  in  the  form 
X(t)  = -|~%]  5U,  (3.23) 


is  kinematically  similar  to 

■ [Hi] 


by  the  transformation 

• [Hi]  ^ 


(3.24) 


(3.25) 


where  X is  a solution  to  Eq  (3.13). 
Proof.  Let 


(3.26) 


Differentiating  Eq  (3.25)  and  substitution  into  Eq  (3.23)  yields 

£+C  [Hi]  ™ - 

Multiplication  of  Eq  (3.27)  by  P 
X as  a solution  to  Eq  (3.13). 

3.3  Reducibility  of  Time  Dependent  Forms  of  Eq  (1.1) 

Theorem  3.5.  The  differential  system  of  Eq  (1.1)  in  the  form 


[Hi]  <3-27> 

reduces  Eq  (3.23)  to  Eq  (3.24)  with 


(3.28) 


where  C(t)  is  an  n x n matrix  whose  elements  are  continuous  real  functions 
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of  a real  variable  and  A,  B are  as  defined  in  theorem  3.3,  is  kinemati- 
cally similar  to 


(3.29) 


by  the  transformation 


X(t)  - [§  ly]  Y(t) 


(3.30) 


where  X(t)  is  a solution  to 


X(t)  *=  AX(t)  - X(t)B  - C(t) 


Proof.  Let 


X(t)  - P Y(t) 


where 


[fPF]  , - [fH^] 


(3.31) 


(3.3) 


(3.32) 


Differentiating  Eq  (3.3)  yields 


X(t)  = 


Y(t)  + 


0 |X(t) 
0 0 


(3.33) 


Substituting  Eq  (3.3)  into  Eq  (3.28)  and  equating  with  Eq  (3.33)  yields 


A |-C(t) 
0 B 


Y(t)  = 


Y(t)  + ° |X~)J  Y(t) 
(3.34) 


Multiplication  and  then  the  addition  of  factors  of  Y(t)  yields 
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AX(t)  - C(t)  - X(t) 
B 


Y(t) 


(3.35) 


J L”  1 A J 


Multiplication  of  Eq  (3.35)  on  the  left  by  P ^(t)  yields 


Y(t) 


AX(t)-X(t)B-C(t)-X( 
B 


X(t)j 


Y(t) 


(3.36) 


which  is  the  desired  result  with  X(t)  a solution  to  Eq  (3.31). 

The  advantage  of  using  this  reduction  technique  can  probably  be 
best  seen  with  an  example  that  will  first  solve  Eq  (1.1)  by  a classical 
method  and  then  by  theorem  3.5. 

Example  3.3. 


(3.37) 


Eq  (3.37)  may  be  written  in  the  following  form: 


x1(t)  = x1(t)  + x2(t)  sin 

t 

(3.38) 

x2(t)  = x2(t) 

(3.39) 

Solving  for  x2(t)  yields 

C x (t)  C 

t { *,<T)  dT  ' t' 

dx 

(3.40) 

Integrating  Eq  (3.40)  and  solving 

for  x^(t)  yields 

x2(t)  = x2(tQ)et  to  = et 

(3.41) 

Substituting  this  result  in  Eq  (3.38)  yields 
x1(t)-x1(t)  = x2(tQ)et  to  sin  t 


(3. 42) 
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using  e~  ps  pn  integrpting  fpctor  in  Eq  (?.42)  yeilds 

C -T  C -t 

/ e x,(T)d  = / x9(t  )e  0 sin  xdt 

t 1 t 2 ° 

o o 


(3.43) 


Integrating  both  sides  and  solving  for  x^(t)  yields 


x^(t)  = 2et-etcos  t 


(3.44) 


XU)  - [2e  ' l cos 

e J 


(3.45) 


The  same  example  using  theorem  3.5  would  be  solved  as  follows:  Let 


A=l;  B=1 ; -C(t)  = sin  t 


(3.46) 


Solving  Eq  (3.31)  for  X(t)  yields 


X(t)  = sin  t 


(3.47) 


Eq  (3.47)  implies  that 


X(t)  = -cos  t 


(3.48) 


except  for  a constant,  and  thus  by  theorem  3.5  Eq  (3.37)  will  reduce  to 


[o  ;] 


(3.49) 


with  Y(t  ) = P (t  ) X(t  ).  Thus 

O O O 


(Wl  [21 

L1. 


(3.50) 


A solution  to  Eq  (3.49)  is 
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yx(t) 


Ly2(t)J 

Then  using  Eq  (3.3)  yields 


"w‘1  . r2eti 

_y2(Co)etJ  LeCJ 


(3.51) 


■ [r^p]  [2i]  - [2et  ■ :! cos '] 


(3.52) 


Although  the  above  example  illustrates  the  two  methods  of  solving 
Eq  (3.28),  the  real  advantage  in  using  theorem  3.5  is  that  it  will  allow 
one  to  solve  higher  order  matrix  equations  when  the  time  dependent  ele- 
ments can  be  arranged  in  one  block.  The  next  definition  will  be  necessary 
in  order  to  show  a technique  that  can  be  used  to  solve  a higher  order 
matrix  equation  in  the  form  of  Eq  (3.28). 


Definition  3.1.  The  tensor  product  ( © ) of  two  matrices  A 


m x m 


and  B is  defined  as 

n x n 


ra..B  a10B  ...  a,  B 
11  12  lm 


A ® B = 


a21B  a22B  •••  a2mB 


La  ,Ba  ,B  ...  a B 
mi  mz  mm 


mn  x mn 


(3.53) 


where 


’ an 

a12 

alm" 

-bll 

b12 

i 

c 

H 

.jo 

A = 

a21 

a^2  • • • 

a2m 

B = 

b21 

b22 

C 1 

CM 

1 

• 1 

• 1 

' aml 

3l  /%  • • • 

m2 

a 

ram  J 

> 

-bnl 

b ■> 
n2 

J 

nn 

(3.54) 
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With  the  above  definition  one  can  readily  show  that  by  using  the  tensor 
product  of  matrices  a one-to-one  (1—1)  correspondence  of  Eqs  (3.13) 
and  (3.31)  can  be  written  in  the  following  form: 


Eq  (3.13)  would  transform  into 


( (A  ® I)  - (I  ® B;  ] X = C = 


Eq  (3.31)  would  transform  into 


fll 

‘in 

?21 

I.*3  nn 


(3.55) 


X(t)  «[  (A  ® 1)  - (I  ® B1)  ] X(t)  - C(t) 


(3.56) 


The  next  two  examples  will  show  how  Eqs  (3.55)  and  (3.56)  would  be 
implemented. 

Example  3.4.  Trying  to  solve  Eq  (3.13)  of  example  3.2  by 
Eq  (3.55)  would  yield 


c h :h:  :h:  j 


'xll 

T 

X12 

3 

X21 

1 

^22" 

.2. 

Multiplication  of  Eq  (3.57)  yields 


10  0 0 
12  0 0 
10  0 0 
0 111 


x. . 

1 

*— 1 

1 

11 

X12 

= 

3 

X21 

1 

_x22  . 

_2_ 

A solution  to  Eq  (3.58)  by  inspection  can  be  readily  seen  as 


(3.57) 


(3.58) 


34 


- 

■ “ 

X11 

1 

X12 

1 

X21 

1 

_x22_ 

0 

lr-1 


C:] 


(3.59) 


The  next  example  will  use  Eq  (3.56)  in  the  process  of  applying 
theorem  3.5. 

Example  3.5.  Solve 


X(t) 


'l 

0 

0 

2 

-t 

e 

t 

-t  -| 
e 

cos  t 

X(t)  ; 

X(0)  = 

ITT 7 

1 

0 

0 

2 

0 

2 

,0 

0 

0 

2 

.2. 

Let 


c a . ■ - to  a , ■ [f  c J 


Solving  Eq  (3.31)  by  use  of  Eq  (3.56)  yields 

-t 


X(t)  = 


-1  0 
0 -1 


0 0 
0 0 


0 0 
0 0 


0 0 
0 0 


X(t)  + 


-t 


t 

cos  t 


A solution  to  Eq  (3.62)  will  then  be  of  the  form 


(3.60) 


(3.61) 


(3.62) 


A(t-t  ) 


X(t)  - e '<>'  7 + ) .i(t-s) 


C(s)  ds 


(3.63) 


where 


a 

A 


-1 

0 

0 

o' 

1 

— — o 

e 

0 

-1 

0 

0 

k-xog. 

1 

C(s)  = 

-s 

e 

0 

0 

0 

0 

2 

s 

_ 0 

0 

0 

0 

1 9 

2_ 

9 

cos 

(3.64) 
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At 

Using  the  methods  of  Chapter  II  to  compute  e , Eq  (3.63)  produces 


X(t)  = 


e C + te  C 
e-C  + te-t 


+ t2/2 
+ sin  t 


(3.65) 


With  X(t)  as  a solution  to  Eq  (3.56),  Eq  (3.60)  reduces  to 


Y(t)  = 


with  Y (t  ) = P-1(t  ) X(t  ).  Thus 
o o o 


1 0 

0 0 

0 2 

0 0 

0 0 

2 0 

0 0 

0 2 

Y (t) 


Y(to) 


0 

1 


0 

0 


-1  -1 

-2  -2 


1 

0 


-3 

-7 

2 


Thus  a solution  to  Eq  (3.66)  is 


Y (t)  = 


-3e 

-7e 

2e 

2e 


2t 

2t 

2t 


(3.66) 


(3.67) 


(3.68) 


Therefore  from  the  substitution  X = P Y,  the  solution  to  the  original 
Eq  (3.60)  is  given  by 


1 

0 

(e  C+te  C) 

(e  C+te  C) 

X(t)  = 

0 

1 

(2  + t2/2) 

(2  + sint) 

0 

0 

1 

0 

0 

0 

0 

1 

-3e 

7 2t 

-7e 

2e2t 

2e2t 


(3.69) 
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or 


'"et  + 4 t eC 

2t  . 2 2t  . - 2t 

e + t e + 2e  sin  t 


S(t) 


2e 


2t 


2e 


2t 


(3.70) 


This  example  demonstrates  a method  that  could  be  used  to  solve  Eq  (3.28). 
At 

The  evaluation  of  e by  the  method  discussed  in  Chapter  II  plays  an 
important  role  in  this  technique. 

3.4  Iterative  Methods  Used  to  Solve  Eg  (3.13) 

This  section  will  examine  possible  iterative  techniques  that 
could  be  used  to  solve  the  Liapunov  equation.  Continuing  research  on 
the  solution  of  this  equation  has  produced  a varied  assortment  of  algo- 
rithms that  could  be  used  to  obtain  a solution  depending  upon  the  proper- 
ties of  the  matrices  A and  B of  Eq  (3.13).  Meirovitch  (Ref  12),  La 
Salle  (Ref  13),  and  Ma  (Ref  14)  give  definite  uses  to  the  solution  of 
Eq  (3.13)  along  with  possible  solution  techniques.  Kitagawa  (Ref  15) 
suggests  an  iterative  technique  for  the  solution  of  an  equation  in  the 
form 

X = F X FT  + S (3.71) 

where  X,  F,  and  S are  n x n matrices.  According  to  Kitagawa  his 
algorithm  works  well  in  terms  of  both  time  and  precision  for  high  dimen- 
sional problems,  and  the  method  is  also  applicable  to  the  Liapunov 
equation  in  the  form 

AT  X + XA  + C = 0 (3.72) 
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The  iterative  technique  in  this  section  was  suggested  by  Lancaster 
(Ref  16)  and  provided  a solution  to  Eq  (3.13)  when  either  or  both  A and 


B were  singular. 

Theorem  3.6.  If 

f (.z)  = (z  + a)  (z  - a)  * (3.73) 

where  a ^ 0 and  a is  real,  and  define 

f (A)  = (aI-A)-1(aI+A)  = U (3.74) 

f (B)  = (al+B) (al-B)-1  * V (3.75) 

then  a solution  of 

X'OXV-f  (U-I)C(V-I)  (3.76) 


is  also  a solution  to  Eq  (3.13). 

Proof.  Multiplication  of  Eq  (3.13)  by  2a  yields 

2aAX  + 2aXB  = 2aC  (3.77) 

2 

Adding  and  subtracting  a X and  A X B to  Eq  (3.77)  yields 

a2X+aXB+aAX+AXB-a2X+aXB+aAX-AXB=2aC  (3.78) 

Factorization  of  Eq  (3.78)  yields 

(aI+A)X(aI+B)-(aI-A)X(aI-B)  = 2aC  (3.79) 

By  choosing  a so  that  (al-A)  ^ and  (al-B)  exist  and  multiplying 
Eq  (3.79)  by  these  factors  on  the  left  and  right  gives 


1 


1 


(aI-A)_1(aI+A)X(aI+B) (al-B)-1  - (aI-A)*1(aI-A)X(aI-B) (al-B)-1  = 

2a(aI-A)C(aI-B)-1  (3.80) 

By  substituting  Eqs  (3.74)  and  (3.75)  into  Eq  (3.80)  yields 

UXV  - X = 2a(aI-A)-1C(al-B)-1  (3.81) 

By  assumption  the  following  equations  are  valid 

f(z)  = il+ii  - f(z)  . x = = 

(z-a)  z-a 

2a(z-a)-1 — f (A)  - I = 2a(A-aI) 

Similarly 


(U-I)  = (A-aI)”1  (3.82) 


(V-I)  = (B-al)-1 

Thus  Eq  (3.81)  transforms  to 

X = UXV  - ^ (U-I)C(V-I) 


where 


(A-al)-1  = -(al-A)-1 
Using  the  recursive  algorithm  yields 

Xn+1  = U Xn  V “ (U_I)  C (V-I)’  (n=1>2»3---> 


(3.83) 


(3.76) 


(3.84) 


(3.85) 


Lancaster  did  not  support  this  iteration  process  with  experimental 
data.  Table  I shows  the  results  of  some  examples  that  were  used  to 
determine  the  usefulness  of  this  algorithm.  The  iteration  process  con- 
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tinued  until  the  answer  was  accurate  to  a specified  number  of  signifi- 
cant digits.  These  examples  were  executed  on  a CDC  6600,  which  is  a 
60  bit  machine,  using  single  precision  arithmetic.  The  most  time  consuming 
part  of  this  algorithm  is  in  determining  the  a parameter  so  that  the 
inverse  of  Eqs  (3.74)  and  (3.75)  exists  for  a converging  sequence 
generated  by  Eq  (.3.85).  The  technique  found  to  work  best  was  as  follows: 

1)  Choose  an  a equal  to  the  minimum  of  the  elements  of  the 
matrices  A,  B,  and  C. 

2)  Assure  that  the  inverse  of  Eqs  (3.74)  and  (3.75)  exists. 

3)  Determine  if  the  sequence  generated  by  Eq  (3.85)  is  a 
converging  sequence. 

4)  If  the  sequence  is  a converging  one, then  iterate  until 
less  than  a specified  tolerance. 

5)  If  the  sequence  diverges  increment  a by  a specified  amount 
and  begin  at  step  2. 

The  specified  amount  suggested  in  step  5 could  be  determined  by  the  . 
elements  of  the  matrices  A,  B,  and  C.  For  example,  if  the  range  of 
the  elements  was  from  .1  to  10  then  an  incremental  value  of  ±1  was 
adequate.  On  the  other  hand,  if  the  range  of  the  values  of  the  elements 
was  from  1 to  100  then  an  incremental  value  of  ±i0  was  adequate. 

Different  solutions  could  be  obtained  by  varying  the  initial  starting 
value  of  XQ.  Also,  with  an  initial  starting  value  of  Xq  = 0 it  was 
usually  possible  to  obtain  more  than  one  solution  by  simply  varying  the 
a parameter.  This  technique  seemed  to  work  well  for  n x n matrices 
where  n < 5.  For  n £ 5 the  technique  is  not  as  consistent  and  a good 
initial  starting  value  is  needed  in  order  to  assure  convergence  to  a 
solution.  For  n > 10  the  iteration  process  breaks  down  due  to  the  amount 
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of  round-off  error  introduced.  Most  10  x 10  matrices  will  not  converge 
even  with  a good  initial  starting  value.  Example  3 of  Table  I is 
for  a 10  x 10  matrix  equation  that  could  not  achieve  better  than  2 


significant  digits  of  accuracy  in  spite  of  a good  initial  starting 
value. 

Another  technique  suggested  by  Lancaster  (Ref  16)  determines  U and 
V as  in  Eqs  (3.74)  and  (3.75)  but  employs  the  following  theorem. 

Theorem  3.7.  An  equation  of  the  form 

X - U X V = Cj  (3.86) 

where 

ci  ■ h (U  " I)  C (V  ' I)  (3,87) 

has  a unique  solution  X with  series  representation 

oo 

X = l uj-1  C.  VJ-1  (3.88) 

j-1 

if  and  only  if  | |u| | • | |v| |<1,  where  | |u| | and  | |v| | are  the  norm  of  the 
matrices  U and  V such  that  the  norm  of  the  identity  matrix  is  1. 

■ 

f 

This  algorithm  proved  difficult  to  implement  since  the  condition 
| | U | |* | | V | | <1  was  not  easily  met  using  a variety  of  norms. 

The  next  theorem  suggests  another  algorithm  that  could  be  used  to 
solve  Eq  (3.13).  The  notation  (*)  means  conjugate  transpose  and  (+) 
means  the  generalized  inverse  of  a matrix. 

Theorem  3.8.  The  series  generated  by 

^ GO 

X * [ l (I-a(I-G)*(I-G))K  a(I-G)*]  (^C)  (3.89) 

k=0 
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4 


where 


0 < a < 


is  a solution  to  Eq  (3.13)  if  and  only  if 


AC  AO 

oTb  is  similar  to  OlB 


(3.90) 


(I-G) (I-G)  (-C)  = -C 


(3.91) 


where 


G - [I  + (A  ® I)  + (I  ® B1)] 


(3.92) 


Proof.  Addition  and  subtraction  of  X to  Eq  (3.13)  yields 


X+AXI+IXB-X=C 


(3.93) 


Using  the  tensor  product  of  definition  3.1  for  factorization  yields 


-X  - (A  ® I)X  - (I  ® B )X  + X = -C 


X - [I  + (A  ® I)  + (I  ® B; ]X  = -C 


(3.94) 


(3.95) 


where  X and  C are  a one-to-one  correspondence  of  the  components  of  the 
matrices  X and  C in  vector  form.  By  substituting  Eq  (3.92)  into 


Eq  (3.95)  yields 


X - GX  « -C 


(l-G)X  » -C 


(3.96) 


(3.97) 
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(3.98) 


X - (I-G)+  (-C) 


Finally,  by  substituting  the  series  representation  for  the  generalized 
inverse  as  per  Ben-Israel  and  Greville  (Ref  17)  one  will  obtain 
Eq  (3.89). 

The  following  example  demonstrates  theorem  3.8. 

Example  3.5.  Solve  Eq  (3.13)  where 


« 1-1  C”1 

‘J  * L°  "2J  • L1 


il 


(3.99) 


Using  Eq  (3.92)  to  generate  C yields 


0 

0 

-1 

0 


0 

0 

0 

-1 


-1 

0 

0 

0 


(3.100) 


Thus 


(I-G)  = 


0 

1 

0 

1 


1 

0 

1 

0 


(3.101) 


Using  the  max  column  norm  for  Eq  (3.101)  and  using  Eq  (3.90)  yields 


0 < a < | 


(3.102) 


Let  a - and  substituting  into  Eq  (3.89)  produces  the  series: 


For  k=0 


-1/2 

-1/2 

-1/2 

-1/2 


(3.103) 
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For  k*l 


X 


(3.104) 


Thus  for  the  2nd  computation 


-1/2 

0 

-1/2 

_1_ 

0 

-1/2 

T 

0 

-1/2 

0 

_ - 

(3.105) 


For  k=3,4,...  the  results  are  the  same  as  for  Eq  (3.105)  and  therefore 


a solution  to  Eq  (3.13)  is 


X 


[-1/2  -1/2*1 
L-l/2  -1/2 J 


(3.106) 


This  method  of  solving  Eq  (3.13)  has  the  advantage  that  one  does  not 
have  to  compute  an  inverse  directly  before  attempting  to  find  a 
solution.  A big  disadvantage,  however,  is  that  it  requires  a large 
amount  of  memory  to  execute  on  a computer  for  n x n matrices  with 
a > 10.- 
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Matrix  A 


Matrix  B = 


Matrix  C = 


Solution  to  Eq  (3.13)  by  Iterative  Method 
Example  1. 

c a 
b a 
q a 


Matrix  X = 
o 


Parameter  a = -2 

Number  of  Iterations  to  3 significant  digits  = 5 


Matrix  Answer  X = 


Matrix  A = 

'5 

7 

9 


[i.oo 

i.ool 

[.25 

•75  j 

Example  2. 

10 

15 

14 

21 

18 

27 

k 


45 


1 


Matrix  B « 


Matrix  B 


1 

2 

4 

8 

16 

32 

64 

128 

256 

100 

0 

1 

2 

4 

8 

16 

32 

64 

128 

50 

0 

3 

6 

12 

24 

48 

96 

192 

175 

120 

100 

256 

128 

64 

32 

16 

8 

4 

2 

1 

50 

206 

73 

14 

21 

32 

77 

125 

85 

95 

7 

14 

21 

28 

35 

42 

49 

56 

63 

70 

9 

18 

27 

36 

45 

54 

63 

72 

81 

90 

10 

20 

30 

40 

50 

60 

70 

80 

90 

100 

0 

10 

20 

30 

60 

70 

80 

100 

110 

150 

1 

8 

36 

57 

28 

56 

150 

90 

24 

71 

Matrix 

C = 

210 

420 

665 

980 

1435 

2170 

3465 

588 

10535 

525 

210 

455 

700 

980 

1330 

1820 

2590 

3920 

6370 

3850 

140 

385 

630 

980 

1540 

2520 

4340 

7840 

7385 

5600 

3745 

9450 

5215 

3220 

2345 

2030 

1995 

2100 

2275 

2485 

1855 

7420 

2870 

910 

1260 

1750 

3430 

5215 

3920 

4375 

525 

1050 

1575 

2100 

2625 

3150 

3675 

4200 

4725 

5250 

385 

770 

1155 

1540 

1925 

2310 

2695 

3080 

3465 

3850 

665 

1330 

1995 

2660 

3325 

3990 

4655 

5320 

5985 

6650 

35 

420 

805 

1190 

2275 

2660 

3045 

3780 

4165 

5600 

[_  385 

980 

2310 

3395 

2730 

4060 

7700 

5950 

3990 

5985 

Matrix 

X 

o 

= 

34 

0 

0 

0 

0 

0 

0 

0 0 0 

0 

34 

0 

0 

0 

0 

0 

0 0 0 

0 

0 

34 

0 

0 

0 

0 

0 0 0 

0 

0 

0 

34 

0 

0 

0 

0 0 0 

0 

0 

0 

0 

34 

0 

0 

0 ( 

3 0 

0 

0 

0 

0 

0 

34 

0 

0 0 0 

0 

0 

0 

0 

0 

0 

34 

0 ( 

3 0 

0 

0 

0 

0 

0 

0 

0 

34  0 0 

0 

0 

0 

0 

0 

0 

0 

0 34  0 

0 

0 

0 

0 

0 

0 

0 

0 0 34 

Parameter 

a = -' 

300 

Number  of  Iterations  to  2 significant  digits  = 4 
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IV.  Reducibillty  Using  Solutions  of  the  Riccati  Equation 


This  chapter  will  extend  the  ideas  of  Chapter  III  to  more  general 
block  forms  of  the  equation 

X(t)  = A(t)  X(t)  (1.1) 

The  first  section  of  this  chapter  will  apply  transformations  to  Eq 
(1.1)  for  the  constant  case  by  solving  the  Riccati  equation  associated 
with  Eq  (1.1).  In  the  second  section  the  discussion  will  cover  specific 
time  dependent  cases.,  and  the  final  section  will  present  some  methods 
that  could  be  used  to  solve  the  Riccati  equation. 

4.1  Reducibility  of  Eq  (1.1)  for  the  Constant  Case 

Theorem  4.1,  The  differential  system  of  Eq  (1.1)  in  the  form 


«*>  -KH]  x<t>  <*•» 

where  A,B,C  and  D are  n x n matrices  belonging  to  the  field  of  real 
numbers,  is  kinematically  similar  to 

- [T  I-dx-bI  ^ <4-2) 

by  the  transformation 

X(7)  = [fjr|]  YU)  (4.3) 

where  X is  a solution  to  the  Riccati  matrix  equation 

AX  + XB  + C + XDX  = 0 (4.4) 

Proof.  Let 

X(t)  = P Y(t)  (4.5) 
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where 


— 1 

p-1  = 

-X  J , 

l n 

io  J 

Eq  (4.5)  yields 

Oil' 
I l-x . 

] Y(t) 

Eq  (4.5)  into  Eq  (4.1)  and  i 

1 l 

1 — 1 

il  *<'>  * [ 

I l-x 

f Eq  (4.8)  on  the  left  by  P 

‘ XD+A 

I -AX-XB-C-XDX  "1 

1 Y(t) 

D 

< -DX-B  J 

-1 


which  reduces  to  Eq  (4.2)  with  X a solution  to  Eq  (4.4). 
Example  4.1.  Reduce  the  differential  system 

X(t) 


Let 


A solution  to  Eq  (4.4)  is 


(4.6) 

(4.7) 

.7)  yields 

(4.8) 

(4.9) 

(4.10) 

(4.11) 

(4.12) 
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Thus  using  Eq  (4.5)  with 


(4.13) 


the  substitution  into  Eq  (4.10)  will  by  theorem  4.1  reduce  the  system  to 


(4.14) 


The  following  theorems  show  the  flexibility  of  block  reducing  Eq  (1.1) 
by  using  Eq  (4.4).  Depending  upon  the  original  differential  system, 
a proper  choice  of  A,B,C,  and  D as  sub-matrices  may  allow  the  solution 
X ir.  Eq  (4.4)  to  be  more  easily  solved.  The  other  alternative  is  to 
vary  the  P matrix  in  Eq  (4.5)  so  that  the  resulting  reduced  differential 
system  is  in  a simpler  form. 

Theorem  4.2.  The  differential  system  of  Eq  (1.1)  in  the  form 


(4.15) 


where  A,B,C,D  are  n x n matrices  belonging  to  the  field  of  real  numbers, 
is  kinematically  similar  to 


Y(t)  = [ 


-DX-B 

-D  -j 

0 

A+XD  J 

Y(t) 


(4.16) 


by  the  transformation 

«3\-  [Hi]  iw 

where  X is  a solution  to  Eq  (4.4). 


(4.17) 
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Proof . As  in  theorem  4.1,  let 


P - 


then  by  differentiating  Eq  (4.17)  and  substitution  into  Eq  (4 


tfrf]  M - [f+§] 


Multiplication  of  Eq  (4.19)  on  the  left  by  P ^ reduces  to  Eq 
with  X a solution  to  Eq  (4.4). 

Theorem  4.3.  The  differential  system  in  the  form  of 
is  kinematically  similar  to 


Y(t) 


A+XD 

0 

• -D 

-DX-B 

Y(t) 


by  the  transformation 


■ [Hr] ?(t) 


where  X is  a solution  to  Eq  (4.4). 
Proof.  Let 


then  by  differentiating  Eq  (4.20)  and  substitution  into  Eq  (4 

trfs]  tm]  - [fHf] 

Multiplication  of  Eq  (4.22)  on  the  left  by  P 1 reduces  to  Eq 
X a solution  to  Eq  (4.4). 


(4.18) 

.15)  yields 

(4.19) 

(4.16) 

Eq  (4.15) 

(4.2) 

(4.20) 

(4.21) 

.15)  yields 

(4.22) 

(4.2)  with 
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Example  4.2.  Reduce  by  theorem  4.3  the  differential  system 


X(t) 


3 0 

1 -2 


4 

4 


-1 

0 


X(t) 


(4.23) 


Let 


[i  .a.  [•;  i],  =-[:  a . »-[s  n 


(4.24) 


A solution  to  Eq  (4.4)  is 


■[-2  i] 


(4.25) 


Using  Eq  (4.21) 


0 -1 
-4  -2 


1 

0 


„-l 


0 

1 


0 

0 


0 

4 


(4.26) 


Eq  (U.23)  by  theorem  4.3  reduces  to 


Y(t) 


0 

0 


-1  -1 
-4  0 


Y(t) 


(4.27) 


The  following  theorems  expand  the  ideas  already  presented  to  other 
forms  of  Eq  (1.1)  by  using  other  partitioning  schemes.  As  in  the 
other  cases  there  is  flexibility  in  what  one  calls  the  matrices  A,B,C, 
and  D and  also  in  the  choice  of  the  P matrix. 

Theorem  4.4.  The  differential  system  of  Eq  (1.1)  in  the  form 
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’-A 

0 

C 

0 

I 

0 

-D 

0 

B 

(4.28) 


is  kinematically  similar  to 


Y(t)  = 


-A-XD 

0 

0 

y 

0 

I 

0 

Y(t) 

(4.29) 

-D 

0 

DX+B 

where  A,B,C,  and  D are  n x n matrices  belonging  to  the  field  of  real 
numbers  by  the  transformation 


9 

I 

0 

-x"j  7 

Kt)  = 

0 

I 

0 Y(t) 

(4.30) 

0 

0 

I 

where  X is  a solution  to  Eq  (4.4), 
Proof.  Let 


I 

0 

-X 

I 

0 

X 

P = 

0 

I 

0 

p 1 = 

0 

I 

0 

0 

0 

I 

9 

0 

m 

0 

I 

(4.31) 


then  by  differentiating  Eq  (4.30)  and  substitution  into  Eq  (4.28)  yields 


I 

0 

-X 

y 

I 

0 

-X 

7 

0 

I 

0 

Y(t)  = 

0 

I 

0 

Y(t)  (4.32) 

0 

0 

I 

0 

0 

I 

Multiplication  of  Eq  (4.32)  on  the  left  by  P reduces  Eq  (4.29)  with 
X a solution  to  Eq  (4.4). 


Theorem 

4.5. 

The  < 

’b 

0 

-D 

X(t)  = 

0 

I 

0 

C 

0 

-A 

(4.33) 


is  kinematically  similar  to 


Y(t)  = 


by  the  transformation 


B+DX 

0 

-D 

0 

I 

0 

Y(t) 

(4.34) 

0 

0 

-XD-A 

I 

0 

0"]  y 

0 

I 

0 Y(t) 

(4.35) 

-X 

0 

I 

where  X is  a solution  to  Eq  (4.4). 

Proof.  The  proof  is  similar  to  theorem  4.4  with 


I 

0 

0 

I 

0 

0 

P - 

0 

I 

0 

p 1 = 

0 

I 

0 

(4.36) 

-X 

0 

I 

* 

X 

0 

I 

Example  4.3.  Reduce  the  differential  system 


X(t)  = 


0 0 

_2 0_ 

0 ] 

0 0 
0 0 


0 0 


(4.37) 


0 -2 


Using  theorem  4.4,  let 

n.  -[o-a],  c.[j  a. 

-D  - [°  J]  (4.38) 

A solution  to  Eq  (4.4)  is 

x * [-4  I2]  (4.39) 

and  thus  by  substitution  for  X into  Eq  (4.29),  Eq  (4.37)  reduces  to 
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Y(t) 


0 -1 

0 0 

0 0 

1 0 

0 0 

0 0 

0 0 

1 0 

0 0 

0 0 

0 1 

0 0 

0 0 

0 0 

0 1 

0 1 

0 0 

4 0 

Y(t) 


(4.40) 


4.2  Reducibility  of  Time  Dependent  Forms  of  Eg  (1.1) 

This  section  will  show  how  the  same  basic  techniques  already 
used  for  constant  matrices  is  also  applicable  to  matrices  with  variable 
elements.  Although  various  choices  of  the  sub-matrices  A,B,C,  and  D 
is  possible  by  choosing  the  proper  transformation  P(t)  matrix,  the 
form  of  Eq  (1.1)  to  be  examined  in  this  section  is 


X(t) 


-B(t) 

D(t) 

-C(t) 

A(t)_ 

X(t) 


(4.41) 


where  A,B,C,  and  D are  n x n matrices  with  polynomial  elements  and  real 
coefficients  such  that 


A(t)  = An  + A,t  + A,t  + ...  + A tP 
Oi^  p 


B(t)  = B + B..  t + B.t  + ...  + B tq 
• U 1 l.  q 


c(t)  = c0  + cxt  + c2t2  + ...  + crtr 


D(t)  = Dq  + Dxt  + D2t2  + ...  + DstS 


(4.42) 

(4.43) 

(4.44) 

(4.45) 


Theorem  4.6.  Let  A,B,C,  and  D be  n x n matrices  with  polynomial 
elements  as  defined  in  Eqs  (4.42),  (4.43),  (4.44),  and  (4.45)  and  let 


X(t  ) = X . Then  a power  series  solution  to 
o o 


X(t)  = A(t)X(t)  + X(t)B(t)  + C(t)  + X(t)D(t)X(t) 


(4.46) 
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is 


X(t)  - X + X t + X-t  + ...  X tn  (4.47) 

o l i.  n 

Proof.  Substitution  of  Eq  (4.47)  into  Eq  (4.46)  yields 

[A  +A1t+A„t2+...+A  tP][X  +X,t+. ..+X  t")  + [X  +X.t+...+X  tn] 
o 1 2 p ol  n 1 o 1 n 

[B  +B1t+...+B  tq]  + [C  +C,t+. ..+C  tr]  + [X  +X.t+. ..+X  tn] 
o 1 q o 1 r o 1 n 

[Dq+D  t+. ..+DgtSJ[Xo+X1t+...+Xntn]  = [Xx+2X2t+3X3t2+. ..+nXntn-1] 

(4.48) 

Multiplication  of  Eq  (4.48)  and  matching  powers  of  t yields  the  follow- 
ing set  of  equations 


A X +X  B +C  +X  D X = X, 

OOOOOOOO  1 


A X, +A, X +X  B1+X,B  +C.+X  D X,+X  D,X  +X,D  X = 2X_ 
ol  lo  ol  lo  1 ool  olo  loo  2 


V2^WAWitV.«2+V.WiW!,1»t 


X, D X, +X, D, X +X„D  X - 3X. 
lol  1 1 o 2oo  3 


(4.49) 

(4.50) 


(4.51) 


The  matching  process  would  then  continue  until  X was  matched. 

n 

Theorem  4.6  can  be  extended  under  certain  conditions  to  include  the 
case  where  A(t),  B(t),  C(t)  and  D(t)  have  elements  as  functions  which 
can  be  approximated  by  a polynomial. 

Theorem  4.7.  The  differential  system  of  Eq  (4.41)  is  kine- 
matically similar  to 
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X(t)D(t)+A(t)  | 0 


-D(t)X(t)-B(t) 


(4.52) 


by  the  transformation 


0 I 

1 -X(t) 


(4.53) 


where  X(t)  in  Eq  (4.47)  is  a solution  to  Eq  (4.46). 


Proof.  Let 


0 I 

1 -X(t)  , 


P_1(t) 


X(t)  I 
I 0 


(4.54) 


Differentiating  Eq  (4.53)  yields 


i(t) ' ffey  + Ftfey  Y(t) 


(4.55) 


Equating  Eq  (4.41)  with  Eq  (4.55),  and  substitution  of  Eq  (4.53)  into 
Eq  (4.41)  yields 


f-B(t)  I D(t)~|  f0  I I I 

L-C(t)  I A(t)J  Li  |-X(t) J 


0 | I 

1 -X(t). 


Y(t)  + 


[rhfuy]  Y(t> 


(4.56) 


Multiplication  and  then  the  addition  of  factors  of  Y(t)  yields 


D(t)  -B(t)-D(t)X(t)  1 _ 0 1 Y(t)  (4  57) 

A(t)  -C(t)-A(t)X(t)+fc(t)  YU)  I -X(t)  C } ' 


Multiplication  on  the  left  by  P (t)  yields 


|~X(t)D(t)+A(t) 
[ D(t) 


-X(t)B(t)-X(t)D(t)X(t)-C(t)-A(t)X(t)+X(t) 

-D(t)X(t)-B(t) 


(4.58) 
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which  reduces  to  Eq  (4.52)  with  X(t)  a solution  to  Eq  (4.46). 


The  next  example  will  apply  theorems  4.6  and  4.7. 
Examp le  4.4.  Reduc  e 


1 

-1 

1 

0 

0 

1 

0 

1 

(t2-t-l) 

(2t-3) 

t 

(t+1) 

_(-2t2+2t+l) 

(-3t+4) 

2 

t 

(4.59) 


(4.60) 


D Tl. 

-G ;] 


Choose 


= f1  c = r (~t2+t+D  (-2t+3)*i 

L°  _1J  * L(2t2-2t-l)  (3t-4)  J, 


(4.61) 


(4.62) 


X('o>  - [o  -i]  ‘ Xo  “•«> 

Applying  Eq  (4.49)  for  a solution  to  Eq  (4.46)  yields 

[; ;]  [i [: ;] 


■>■[.:  a 


(4.63) 

(4.64) 
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Similarly  applying  Eqs  (4.50)  and  (4.51)  yields 
*2  • X3  ' [o  o] 


Thus  using  Eq  (4.47)  a solution  to  Eq  (4.46)  is 


(4.65) 


IX  -i] 


(4.66) 


Therefore  by  theorem  4.7  with 


P(t)  = 


0 0 10 

0 0 _0 1 

1 0 -t  -1 

0 1 t 1 


(4.67) 


Eq  (4.59)  reduces  to 


2t  (t+2) 

0 

0 

(2-t) (t-1) 

0 

0 

1 0 

(-t-1) 

-2 

0 1 

t 

2 

(4.68) 


with  Y(tQ)  = P (to)X(to),  so  that 


Kt0)  - 


(4.69) 


The  differential  system  of  Eq  (4.68)  can  now  be  reduced  still  further 
by  solving  Eq  (4.46)  with  D(t)  = 0 and  using  the  following  theorem. 
Theorem  4.8,  The  differential  system  in  the  form 

«*>  ' [l-gctj  ]'A(t>]  (4-70> 

where  A,B,  and  C are  defined  in  Eqs  (4.42),  (4.43),  and  (4.44),  is  kine- 
matically similar  to 
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Y(t) 


r -B(t)  o i 
! 0 I A(t)J 


Y(t) 


with  the  transformation 

I 


X(t) 


-X(t) 


Y(t) 


where  X(t)  is  a solution  Eq  (4.46)  with  D(t)  = 0. 
Proof.  Let 


P(t) 


-X(t)  I 


P_1(t) 


X(t) 


(4.71) 


(4.72) 


(4.73) 


As  in  theorem  4.7,  differentiating  Eq  (4.72),  equating  with  Eq  (4.70) 
with  the  substitution  of  Eq  (4.72)  yields 


-B(t) 

0 

_ -C(t) 

A(t)_ 

-X(t) 


Y(t)  = 


-X(t) 


Y(t)  + 


0 

1 0 

-x(t) 

1 0 

Y(t) 


(4.74) 


Then  multiplication  on  the  left  by  P ^(t)  and  grouping  the  factors  of 
Y(t)  yields 


-B(t) I 0 1 

-A(t)X(t)-X(t)B(t)-C(t)+X(t)  I A(t)J 


Y(t) 


(4.75) 


which  reduces  to  Eq  (4.71)  with  X(t)  a solution  to  Eq  (4.46)  with 
D(t)  = 0. 

One  can  readily  see  that  the  process  could  be  carried  still  further 
by  separating  the  system  of  Eq  (4.71)  into  two  separate  systems  and 
then  alternating  theorems  4.7  and  4.8  until  a total  reduction  was  obtained. 
The  next  section  will  demonstrate  this  technique  for  the  constant  case. 
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This  section  will  examine  solution  techniques  to  solve  the  Riccati 


equation.  As  with  the  Liapunov  equation,  new  methods  to  solve  various 
forms  of  Eq  (4.4)  are  under  constant  study.  The  Riccati  equation  by 
itself  plays  an  important  part  in  the  optimal  estimation  and  filtering 
problem  pursued  by  Kalman  (Ref  18)  and  also  in  orbital  trajectory  opti- 
mization as  investigated  by  Bryson  and  Ho  (Ref  19),  Breakwell  (Ref  20), 
and  McReynolds  (Ref  21).  As  was  mentioned  in  Chapter  I,  a thesis  by 
Leuthauser  (Ref  4)  investigates  ways  to  obtain  solutions  to  Eq  (4.4), 
and  this  section  will  apply  some  of  the  methods  suggested  by  Leuthauser. 

Probably  one  of  the  best  existence  theorems  in  obtaining  a solution 
to  Eq  (4.4)  was  done  by  J.  Jones  (Ref  22),  and  it  has  a direct  applica- 
tion to  the  theorems  of  this  chapter. 

Theorem  4.9.  (Jones  (Ref  22)).  Let  f_^(A)  be  any  polynomial 
of  degree  n>l  in  A with  coefficients  belonging  to  the  field  of  real 
numbers  such  that 

VR>  ■ [fH]  <4-80> 

where  A,B,C,  and  D are  of  order  n x n, 
then  a solution  of 

(X,I)f_^(R)  = (0,0) 

with  U 1 or  M ^ existing,  or  a solution  of 

V*>  [-x>  [o] 

with  M ^ or  N ^ existing,  is  also  a solution  Eq  (4.4). 


(4.81) 


(4.82) 
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I 


Theorem  4.9  when  combined  with  theorem  4.1  can  constitute  one  theorem 
as  was  originally  presented  by  Jones,  Lukes,  and  Dolan  (Ref  23).  The 
next  example  will  use  theorem  4.1  to  reduce  a differential  system  and 
also  theorem  4.9  along  with  the  work  done  by  Leuthauser  to  solve 
Eq  (4.4). 


Example  4.5.  Reduce  the  differential  system 


(4.83) 


Substitution  of  R for  the  matrix  of  Eq  (4.83)  and  finding  the  characteristics 
roots  yields 


| XI-R|  = A4-5A3+9A2-7A+2=0 


A-l, 1,1,2 


(4.84) 


Constructing  the  table  of  characteristic  roots  as  suggested  by  Leuthauser 
yields 


TABLE  II 


Combinations  of  Characteristic  Roots  of  R 


1 

1 

1 

2 


(A-l)2 
(A-l)2 
(A-l)2 
(A-l) (A— 2) 


(A-l)2 
(A-l)2 
(A-l)2 
(A-l) (A— 2) 


(A-l)2 

(A-l)2 

(A-l)2 

(A-l)(A-2) 


2 

(A-l) (A— 2) 
(A-l) (A— 2) 
(A-l) ( A— 2 ) 
(A— 2) 2 


Table  II  shows  that  there  are  only  3 possibilities  to  check  for  a solu- 
tion to  Eq  (4.4):  (A-l)2,  (A-l) (A-2) , and  (A— 2) 2 . Choosing  (A-l)2  as 
a possible  solution 


] 


63 


fx(R)  - (R-I)2  = R2-2R  + I (4.85) 

Substituting  the  matrix  R into  Eq  (4.85)  and  using  Eq  (4.80)  yields 


(4.86) 

Solving  the  equations 


XU  + V = 0 (4.87) 

XM  + N = 0 (4.88) 


from  Eq  (4.81)  yields  from  solving  Eq  (4.87) 


(4.89) 


A check  of  this  solution  with  Eq  (4.88)  confirms  by  theorem  4.9  that  X 
is  a solution  to  Eq  (4.4).  Thus  with  one  valid  solution  X,  by  theorem 
4.1  the  system  of  Eq  (4.83)  reduces  to 


1.34 

.33 

0 

0 

-.32 

.66 

0 

0 

-1 

-4 

3 

2 

5 

-4 

6 

2 

Y(t) 


(4.90) 


In  both  Chapters  III  and  IV,  the  kinematic  similarity  theorems  have 
been  limited  to  the  case  where  A,B,C,  and  D were  n x n matrices.  All 
these  theorems  are,  however,  just  as  valid  for  the  partitioning  of  the 
sub-matrices  A,B,C,  and  D into  rectangular  blocks.  The  only  real  require- 
ment is  that  the  resulting  Liapunov  or  Riccati  equation  is  conformable 
to  obtaining  a solution.  Although  theorem  4.9  as  originally  done  by 
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J.  Jones  was  specifically  for  the  case  where  the  sub-matrices  were  of 
order  n x n,  the  theorem  can  be  extended  to  cover  rectangular  matrices. 


Theorem  4.10.  Let  f^(X)  be  any  polynomial  of  degree  n>l  in 
X with  coefficients  belonging  to  the  field  of  real  numbers  such  that  R, 
fy(R)  as  given  by  Eq  (4.80)  for  A^,  Bpxp,  C^,  DpxM-  Then  a solu- 
tion of  Eq  (4.81)  with  U ^ existing  or  a solution  to  Eq  (4.82)  with 
N l existing  for  X^^^,  is  also  a solution  to  Eq  (4.4). 

Proof.  Similar  to  the  proof  done  by  J.  Jones  (Ref  22)  one 
must  first  show  that  if  a solution  of  Eq  (4.31)  exists  with  U 1 
existing  then  a solution  is  also  a solution  of  Eq  (4.4).  Let  such  a 
solution  X exist.  Now  the  matrices  R,  f (R)  commute  so  we  have  the 
following  identities: 

DV+UB  = BU  - MC;  UD  + MA  = DN  - BM  (4.91) 

AV+VB  = CU  - NC;  CM  + VD  = -NA  + AN 

Making  use  of  Eqs  (4.91)  and  the  above  solution  X yields 

0 = -(XM+N)C 
0 = -XUB-XMC+XUB-NC 
0 =*  -X(UB+MC)+XUB-NC 
0 = -X(BU-DV)+XUB-NC 

0 = -XBU-XDXU+XUB-NC  (4.92) 

0 = -XDXU-AXU-XBU+AXU+XUB-NC 
0 = -XDXU-AXU-XBU-AV-VB-NC 
0 = -XDXU-XBU-AXU-CU 
0 = - ( XDX+AX+XB+C ) U 

and  since  U * exists  X is  a solution  to  Eq  (4.4). 

Next  let  X be  a solution  of  Eq  (4.82)  with  N ^ existing.  Making 
use  of  Eqs  (4.91)  yields 
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f 


0 - -A(V-NX) 

0 - (-CU+NC+VB) + (CM+VD+NA) X 
0 * -CMX+NC+VB+CMX+NXDX+NAX 

0 = -AV+CU+NXnx+NAX  (4.93) 

0 * VB+NC+NXDX+NAX 
0 = NXB+NC+NXDX+NAX 
0 - N (XDX+AX+XB+C) 


Thus  Eqs  (4.93)  imply  that  X is  a solution  to  Eq  (4.4). 

The  next  example  will  demonstrate  how  the  reducibility  theorems 
can  implement  a rectangular  partitioning  scheme,  and  it  will  also  show 
how  the  alternate  use  of  theorems  4.1  and  theorem  3.4  can  reduce  a 
differential  system  like  Eq  (1.1)  to  diagonalized  systems. 
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1 


' 2 -4"|  -C 

r5 

- (2,2),  A = 6 

(4.98) 

A solution  to 

Eq  (3.13)  is 

X = 

[T] 

(4.99) 

Therefore  Eq 

(4.97)  reduces 

to 

Z(t) 

. r 2 -4  oi 

Ut) 

(4.100) 

-5  10 

0 0 6 


Eq  (4.100)  can  now  be  broken  into  two  separate  systems 


(4.101) 


(4.102) 


Applying  theorem  4.1  on  Eq 


(4.101),  let 


-B=2,  -C=-5,  D=-4,  A=1 


A solution  to  Eq  (4.4)  is 


X = 1 


Thus  Eq  (4.101)  reduces  to 


Using  theorem  3.4  again  yields 


(4.103) 


(4.104) 


(4.105) 
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t 


B— 3,  -C=-4,  A=6 


(4.106) 


Solving  Eq  (3.13) 


X * 4/9 


(4.107) 


Therefore  Eq  (4.105)  reduces  to 


S(t)  = 


[1  s] 


(4.108) 


Thus  by  saving  the  solutions  to  Eqs  (4.4)  and  (3.13),  a solution  to 
the  original  differential  system  of  Eq  (4.94)  can  be  found  by  solving 
the  diagonalized  systems  of  Eqs  (4.108)  and  (4.102). 

By  examining  the  matrix  used  in  this  example,  which  was  found  in 
Nerring  (Ref  7:216),  the  classical  reduction  technique  for  reducing 
this  system  would  require  finding  the  generalized  eigenvectors  in  order 
to  reduce  the  matrix  to  a Jordan  form. 

Iterative  Methods  to  Solve  Eq  (4,4),  The  iterative  methods 
investigated  in  this  section  follow  from  the  work  done  by  Leuthauser 
(Ref  4),  and  it  is  an  extension  of  Lancaster's  method  (Ref  16)  applied 
to  Eq  (4.4). 

Theorem  4.11  (Leuthauser  (Ref  4)).  If 


f(z)  = (z+a)(z-a) 


where  a ^ 0 and  a is  real. 

and  if  U and  V are  as  defined  in  theorem  3.6  then  a solution  of 


X =■  UXV  - (U-I)C(V-I)  - jg  (U-I)XDX(V-I) 


is  also  a solution  to  Eq  (4.4), 


(4.109) 


(4.110) 
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The  process  used  to  determine  the  a parameter  was  identical  to  the 


one  discussed  in  Chapter  III  in  finding  a solution  to  the  Liapunov 

equation.  Different  solutions  could  be  obtained  by  varying  the  initial 

starting  value  X . Also,  with  an  initial  starting  value  X = 0 it 
o o 

was  usually  possible  to  obtain  more  than  one  solution  by  simply  varying 

the  a parameter.  Table  III  shows  the  results  of  some  examples  that 

were  executed  using  single  precision  arithmetic  on  the  CDC  6600.  As 

with  the  Liapunov  equation  for  n x n matrices  where  n<5,  solutions  were 

readily  obtained.  However,  for  n>5  a good  initial  starting  value  was 

needed.  For  n£l0  it  is  difficult  to  make  the  iterative  process  converge 

even  with  a good  initial  approximation.  Example  3 of  Table  III  is  for 

10  x 10  matrices  which  could  only  be  made  accurate  to  2 significant 

digits  before  iterating  away  from  a solution.  Of  special  note  in 

Table  III  is  Example  2 which  would  not  yield  a solution  by  the  method 

of  Example  4.5  because  matrices  U,M,  and  N were  all  singular. 

This  iterative  technique  can  also  find  the  square  root  of  a matrix 

by  setting  A=B=0,  D = -I,  and  C as  the  input  matrix.  Table  IV  shows 

the  results  of  a square  root  example.  For  the  square  root  case  once 

an  a parameter  was  found,  -a  would  yield  the  opposite  sign  solution. 

The  initial  starting  value  of  X =0  seemed  to  work  well  for  n£l0. 

o 

Another  iterative  technique  investigated  was  an  adaptation  of  a 
method  suggested  by  Carnhahan  (Ref  24)  for  solutions  of  equations  in 
the  form 

X = F(X)  (4.111) 

By  letting  F(X)  be  Eq  (4.110)  then  similar  to  Carnhahan  a new  iterative 
scheme  in  the  form 
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X - (I-K)X  + K F(X) 


(4.112) 


where  K is  a matrix  parameter  that  speeds  convergence,  could  be  used  to 
solve  Eq  (4.4). 

The  same  examples  of  Table  III  were  also  tested  on  Eq  (4.112)  with 
very  little  difference  in  the  results  except  for  example  2 which  reached 
an  answer  in  20  iterations. 


70 


TABLE  III 

Solutions  to  Eg  (4.4)  by  Iterative  Method 

Example  1. 


Matrix  A = 

Matrix  B = 

Matrix  C = 

Matrix  D = 

Matrix  X = 
o 


■ 

L 


a 


[! 


i 

0 


] 


0 

1 

0 

0 


Parameter  a = 3 

Number  of  Iterations  to  3 significant  digits  = 4 
Matrix  Answer  X = 


’ .414  -.414 

-.414  .414 


Parameter  a = -7 

Number  x>f  Iterations  to  3 significant  digits  = 8 
Matrix  Answer  X = 


-2.41  2.41 

2.41  -2.41 


Example  2. 


Matrix  A = 


Matrix  B = 


4 0 0 

13  0 

112 

-2  -1  o' 

1-4  0 

1 -1  -2 
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Matrix  C ■ 


Matrix  D = 


Matrix  X 

o 


Parameter  a = -4 

Number  of  iterations  to  3 significant  digits  = 31 


Matrix  Answer  = 

-.939 

.939 

-.939 

-.939 

.939 

-.939 

-.939 

.939 

-.939 

Example  3. 

Matrix  A = 


5 

10 

15 

20 

25 

30 

35 

40 

45 

50 

6 

12 

18 

24 

30 

36 

42 

48 

54 

60 

4 

8 

12 

16 

20 

24 

28 

32 

36 

40 

7 

14 

21 

28 

35 

42 

49 

56 

63 

70 

3 

6 

9 

12 

15 

18 

21 

24 

27 

30 

8 

16 

24 

32 

40 

48 

56 

64 

72 

80 

2 

4 

6 

8 

10 

12 

14 

16 

18 

20 

9 

' 18 

27 

36 

45 

54 

63 

72 

81 

90 

1 

2 

3 

4 

5 

6 

7 

8 

9 

10 

10 

20 

30 

40 

50 

60 

70 

80 

90 

100 

Matrix  B = 


1 

2 

4 

8 

16 

32 

64 

128 

256 

100 

0 

1 

2 

4 

8 

16 

32 

64 

128 

50 

0 

3 

6 

12 

24 

48 

96 

192 

175 

120 

iioo 

256 

128 

64 

32 

16 

8 

4 

2 

1 

50 

206 

73 

14 

21 

32 

77 

125 

85 

95 

7 

14 

21 

28 

35 

42 

49 

56 

63 

70 

9 

18 

27 

36 

45 

54 

63 

72 

81 

90 

10 

20 

30 

40 

50 

60 

70 

80 

90 

100 

0 

10 

20 

30 

60 

70 

80 

100 

110 

150 

1 

8 

36 

57 

28 

56 

150 

90 

24 

71 

72 


Matrix  C = 


r- 

714 

- 9380 

-21385 

-28420 

-35315 

-41930  -47985 

-52920 

-55615 

3500 

- 3465 

-10570 

-14000 

-17395 

-43995 

-181930  -70910 

-106330 

-18130 

1750 

- 4760 

-13090 

-58170 

-72520 

-120960 

-107730  -69160 

-41160 

-11715 

4200 

-118755 

-64050 

-43785 

-21280 

-89530 

-34720  -22505 

-181650 

-22225 

35 

-102270 

-66080 

-43680 

-109340 

-53865 

-71750  -82320 

-20510 

-58555 

3325 

-43575 

-50400 

-29050 

-34650 

-46375 

-70350  -20825 

-106050 

-93275 

2450 

-21665 

-24955 

-42945 

-86660 

1925 

2310  -46305 

-143920 

-70035 

3150 

1- 

560 

- 1120 

- 1680 

- 2240 

-2800 

-3360  4655 

- 4480 

- 5040 

3500 

1-61215 

-24080 

-244195  -151935 

-83475 

-70840  -45955 

-118720 

-32585 

5250 

L 

350 

700 

1050 

1400 

1750 

2100  2450 

2800 

3150 

0 

Matrix 

D = 

6 

8 

18 

24  30 

36 

42 

48  54 

60 

3 

9 

12 

15  37 

150 

60 

90  20 

37 

4 

11 

48 

60  100 

90 

60 

40  20 

70 

100 

60 

40 

20  75 

30 

20 

150  20 

0 

85 

60 

38 

90  45 

60 

70 

21  51 

77 

36 

42 

25 

30  40 

60 

20 

90  80 

31 

18 

21 

36 

72  0 

0 

40 

120  60 

200 

1 

2 

3 

4 5 

6 

0 

8 9 

10 

50 

20 

200  125  70 

60 

40 

100  30 

0 

2 

4 

6 

8 20 

40 

90 

15  7 

6_ 

Matrix 

X = 

0 

34 

0 

0 

0 0 

0 

0 

0 0 

0 

0 

34 

0 

0 0 

0 

0 

0 0 

0 

0 

0 

34 

0 0 

0 

0 

0 0 

0 

0 

0 

0 

34  0 

0 

0 

0 0 

0 

0 

0 

0 

0 34 

0 

0 

0 0 

0 

0 

0 

0 

0 0 

34 

0 

0 0 

0 

0 

0 

0 

0 0 

0 

34 

0 0 

0 

0 

0 

0 

0 0 

0 

0 

34  0 

0 

0 

0 

0 

0 0 

0 

0 

0 34 

0 

0 

0 

0 

0 0 

0 

0 

0 0 

0 

Parameter  a = - 

3000 

Number 

of  iterations  to  2 significant  digits 

= 3 
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TABLE  IV 


Eq  (4.110)  to  Solve  Square  Root  of  a Matrix 


Example  1 


Matrix  A = 0 


10x10 


Matrix  B = 0 


10x10 


Matrix  C = 


Matrix  D = -I 


10x10 


Matrix  X = 0,» 

o 10x10 


Parameter  a = -30 


Number  of  iterations  to  3 significant  digits  = 17 
Matrix  Answer  = 


3.16 

.012 

.022 

.031 

.038 

.045 

.051 

.057 

.062 

.067 

.013 

4.47 

.019 

.027 

.034 

.040 

.046 

.051 

.056 

.061 

.011 

.010 

5.48 

.025 

.031 

.037 

.042 

.048 

.053 

.057 

.010 

.018 

.025 

6.32 

.029 

.035 

.040 

.045 

.050 

.054 

.009 

.017 

.023 

.029 

7.07 

.033 

.038 

.043 

.047 

.052 

.009 

.016 

.022 

.028 

.033 

7.74 

.037 

.041 

.046 

.050 

.008 

.015 

.021 

.027 

.032 

.037 

8.37 

.040 

.044 

.048 

.008 

.015 

.020 

.026 

.031 

.035 

.040 

8.94 

.043 

.047 

.008 

.014 

.020 

.025 

.030 

.034 

.039 

.043 

9.49 

.046 

.007 

.013 

.019 

.024 

.029 

.033 

.038 

.041 

. 046 

10.0 

V.  CONCLUSION 


The  kinematic  similarity  theorems  along  with  methods  to  compute  a 
fundamental  solution  to  the  differential  system  of  Eq  (1.1)  should 
prove  beneficial  for  a wide  variety  of  applications.  Close  scrutiny 
of  the  numerous  examples  should  show  the  adaptability  of  these  algorithms 
to  a digital  computer.  In  particular,  the  implementation  of  these 
theorems  for  the  time  dependent  case  should  be  notably  more  efficient 
than  the  time  consuming  iterative  methods  that  yield  questionably 
accurate  solutions.  For  the  constant  cases  of  Eq  (1.1),  the  application 
of  these  theorems  could  greatly  simplify  the  original  differential  system 
so  that  conventional  solution  methods,  which  might  not  have  otherwise 
worked,  could  be  applied  to  the  kinematically  similar  differential  system. 
Finally,  although  further  research  is  needed  on  the  iterative  methods 
used  to  solve  the  Liapunov  and  Riccati  equations,  these  iterative  scnemes 
are  able  to  find  solutions  that  are  not  readily  found  by  other  methods. 
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